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ABSTRACT 


/  ■  '  ’ 

,We  show-,  that  a  uniformly  distributed  random  deviate  between  zero  and  one 
can  be  constructed  by  combining  two  independent  random  variables:  1)  a  random 
exponent  drawn  from  a  geometric  distribution,  and,  2)  a  random  mantissa  drawn 
from  a  uniform  distribution.  Algorithms  for  generating  these  random  variables 
are  developed,  and  a  Pascal  procedure  combining  these  into  a  computationally 
efficient  random  number  generator  for  micro  computers  is  presented.  The 
resulting  generator  is  computationally  efficient  as  itx  1)  uses  no  floating 
point  arithmetic,  and,  2)  uses  on  the  average  only  3.5  random  bytes  to 

f 

construct  a  four  byte  random  deviate.  Finally  we  she*  empirically  that  the 
generator  has  both  an  unusually  long  period  and  excellent  statistical 
properties . 


AMS  (MOS)  Subject  Classifications:  65C10,  68J10,  90-04 

Key  Words:  Random  Number  Generators,  Uniform  Distribution, 
Micro  Computers 
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SIGNIFICANCE  AND  EXPLANATION 


Conventional  generators  of  uniformly  distributed  deviates  between  zero 
and  one  usually  obtain  their  deviates  by  converting  a  potentially  large  random 
integer  to  floating  point  notation  and  scaling  it  down  as  required.  This  is  a 
reasonable  procedure  for  a  computer  with  a  large  words! ze  and  with  floating 
point  arithmetic  implemented  in  hardware.  However  many  micro  and  mini 
computers  use  16  bit  integer  arithmetic  and  have  floating  point  arithmetic 
implemented  in  software.  For  these  computers,  the  conventional  approach  to 
generation  of  uniformly  distributed  deviates  between  zero  and  one  is  quite 
undesirable  as  it  is  both  time  consuming  and  yields  deviates  with  a  short 
period  and  a  low  resolution. 

In  this  paper  we  present  an  alternative  approach.  The  resulting 
generator  is  computationally  efficient  and  has  good  statistical  properties. 

In  addition  it  has  an  unusually  long  period  (we  stopped  measuring  at 
64,000,000)  even  when  driven  by  short  period  (32,768)  congruential  generators, 
and  like  numbers  appear  at  random  intervals.  The  generator  is  based  on  the 
fact  that  a  floating  point  number  is  represented  in  the  computer  in  terms  of 
two  variables  (a  mantissa  and  an  exponent).  These  are  developed  independently 
of  each  other  without  the  use  of  floating  point  arithemtic. 

The  good  performance  characteristics  of  the  proposed  generator  are 
probably  due  to  the  following:  (1)  no  floating  point  arithmetic  is  performed, 
(2)  a  variable  number  of  random  bytes  (between  three  and  six)  is  used  to 
construct  a  number,  (3)  half  the  time  the  process  is  extremely  simple  and 
requires  only  three  bytes,  and,  (4)  on  the  average  about  3.5  random  bytes  are 
required  to  construct  a  four  byte  random  floating  point  number. 


The  responsibility  for  the  wording  and  views  expressed  in  this  descriptive 
summary  lies  with  MRC,  and  not  with  the  author  of  this  report. 
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Conventional  generator  »  of  tin  i  for  ml  v  distributed  de-.-i  a  tee  between 
aero  and  one  usually  obtain  their  deviates  bv  converting  a 
potentially  large  random  integer  to  floating  point  notation  and 
scaling  it  down  as  required.  This  :ls  a  reasonable  procedure  tor  a 
computer  with  a  large  wordsize  and  with  floating  point  arithmetic 
implemented  in  hardware.  However  many  micro  and  mini  computer  s 
use  16  bit  integer  arithmetic  and  have  floating  point  arithmetic 
implemented  in  software.  For  these  computers,  the  conventional 
approach  to  generation  of  uniformly  distributed  deviates  between 
zero  and  one  is  quite  undesirable  as  it  is  both  time  consuming 
and  yields  deviates  with  a  short  period  and  a  low  resolution. 

In  this  paper  we  present  an  alternative  approach.  The  resulting 
generator  is  computationally  efficient  and  has  good  statistical 
properties..  In  addition  it  has  an  unusually  long  period  (we  stop¬ 
ped  measuring  at  64,000,000)  even  when  driven  by  short  period 
(32,768)  congruent! al  generators,  and  like  numbers  appear  at 
random  intervals.  The  generator  is  based  on  the  fact  that  a 
floating  point  number  is  represented  in  the  computer  in  terms  of 
two  variables  (a  mantissa  and  an  exponent  i.  These  are  developed 
independently  of  each  other  without  the  use  of  floating  point 
ar i thmeti c . 

Our  implementation  of  the  algorithm  is  based  on  the  conventions 
for  representation  of  floating  point  numbers  adopted  by  the  AMD 
9311  iOR  INTEL  8231A)  hardware  floating  point  unit  til.  This  is 
the  convention  used  in  most  micro  computer  compilers  such  as  Pas¬ 
cal  /MT+  and  FORTRAN-80.  This  convention  differs  slightly  from 
(but  is  conceptual l v  identical  to)  the  proposed  IEEE  standard 
notation  used  in  many  other  systems  C23 .  The  decision  to  use  the 
de  facto  micro-computer  standard  rather  than  the  IEEE  standard 
reflects  the  present  unavai 1 abi 1 1 t v  of  efficient  floating  point 
random  number  generators  on  micro  computers.  An  implementation  of 
the  algorithm  under  the  IEEE  standard  should  follow  our  examples 
quite  closely. 

The  good  performance  character i s L 1  re  ot  the  proposed  generator 
are  probably  due  to  the  following:  (1)  no  floating  point  arit¬ 
hmetic  is  performed  *2)  a  variable  number  of  random  bytes 
(between  three  and  six)  is  used  to  construct  a  number  .  CV>  hul  f 
the  time  the  process  if  extremely  s  1  mo  1  e  arid  requires  on  l  v  three 
bvtes,  and,  \4>  on  the  average  about  3.3  random  bvtes  are  re¬ 
quired  to  construct  a  tour  bvte  random  floating  go '  nt.  number  . 

1- 1 —  o  &  r  x  1 -1 43  f-  o  1  i  -j  r  m  o  r  a  r  i  c:  >  dm 

All  floating  point  conventions  are  based  on  the  fact  that  anv 
rational  number  X  cun  be  r  «pr  e  son  ted  as  l  he  prodnc  I  ot  a  fractio¬ 
nal  part  and  an  integer  part.  At  the  cos -I  t  t,ol  bei  ug  able  t.o 
accurately  represent  all  rational  numbers.  these  convent  ions  are 
j  mu  1  emen t e< f  ■..< s  fid  1  ows : 
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*  2 


Equation  1 


M 


where  X  =  The  number  to  be  represented 

•  *  an  integer  exponent 

h  =  a  constant  (usually  a  power  of  two). 

m  =  the  mantissa,  m  >=M/2  unless  X  =  0  in  which  case  m  =  0 

The  precision  of  X  is  determined  by  M  while  the  range  of  X  is 

determined  by  the  permissible  range  of  e.  Note  that  the 

fractional  part  <m/M>  of  X  is  always  in  the  range  1/2  to  1 
regardless  of  the  value  of  X  (as  long  as  X  i s  nonzero).  In  Figure 
1  we  present  the  ranges  of  values  taken  on  by  e  and  m/M  for 
different  values  of  X  between  zero  and  one. 

+ - + — + - 1 

•  I  !  •  !  a/H  ! 

♦  - ♦ — ♦ - ♦ 

•1.0  >x>=0.5  !0  10.5-1.0! 

♦ - * — ♦ - ♦ 

•0.5  >x>=0.25  !-l  !  0.5  -  1.0  ! 

•  f - ♦ — ♦ - ♦ 

! 0. 25  >*>=0.125  !-2  0.5  -  1.0  ! 

♦ - ♦ — ♦ - ♦ 

10.125  >*>=0.0425  !-3  !  0.5  -  1.0  ! 

!0.0425>* >=0.03125! -4  !  0.5  -  1.0  ! 

♦ - ♦ — ♦ - ♦ 

Figure  Is  Values  of  •  and  m/M  for  numbers  on  CO-11 

It  is  seen  that  the  exponent  e  defines  the  range  of  a  number 

(i.e.  0.5  >  x>=0.25>,  and  that  the  mantissa  m  defines  the 

specific  value  of  the  number  within  this  range. 

RANDOM  VARIABLES  OIM  CO-1  □ 

From  the  preceding  discussion,  we  see  that  a  random  variable  u  on 
CO-13  can  be  expressed  as  a  function  of  a  random  exponent  e  and  a 
random  mantissa  m  as  follows: 

m  e 

u  =  -  *  2  Equation  2 

M 

where  u  =  random  variable  on  CO-11 

•  =  random  variable  drawn  from  the  distribution: 

-(i  +  1 ) 

Prim  =  i  3-  =  2  i  =  0,1,2,.. 

M  =  a  constant 

m  -  Uniformly  distributed  random  variable  on  C  M/2,M  ). 


Perhaps  the  most  useful  consequence  of  the  fact  that  u  is 
expressed  as  a  function  of  two  independent  random  variables  is 

the  fact  that  it  enables  us  to  propose  an  algorithm  with  a  period 

so  long  that  we  have  been  unable  to  empirically  measure  it. 

As  shown  in  Figure  2,  the  algorithm  starts  out  by  assigning  a 
value  of  zero  to  the  exponent  (defining  a  number  on  CO. 5-1. 03) 
and  a  random  value  to  the  mantissa.  A  check  is  made  to  see  if  the 
resulting  mantissa  has  a  valid  value  (there  is  a  50 7.  probability 
that  this  is  so),  if  it  does,  the  algorithm  stops  as  a  valid 
number  has  been  produced. 

1.  INITIAL  ASSI6NHENTS 

e  -  0  (Correct  value  if  u  >  1/2) 

■  *  U<0,N)  (Uni fora  deviate  on  0-N) 

II.  IS  ADDITIONAL  MAX  REQUIRED? 

If  a  >=  H/2  (  this  happens  half  the  tiee) 

then  goto  step  IV 
else  a  z  a  +  H/2. 

III.  ADJUST  EXPONENT 

A.  Draw  randoe  byte(s)  until  a  nonzero  byte  is  found: 

k  *  RandoeDyte 

uhile  k  =  0  (  this  happens  uith  a  probability  of  1/256  ) 

e  »  e  -  8 
k  -  RandoeDyte. 

B.  Scan  the  byte  for  the  first  nonzero  bit: 

uhile  k  <  120 
k  =  k  I  2 
e  =  e  -  1. 

IV.  RANDON  VARIABLE  IS  u  -  He, el 

Figure  2:  Algorithm  for  generating  u 


If  the  mantissa  is  not  valid  (i.e.  its  first  bit  is  not  one), 
then  M/2  is  added  to  m  to  make  it  valid,  and  a  procedure  for 
generating  a  random  exponent  is  entered.  This  procedure  is  based 
on  the  premise  that  the  number  of  zeroes  preceding  the  first  one 
in  a  random  stream  of  bits  follows  the  geometric  distribution 
with  p  =0.5.  Random  bytes  are  drawn  until  a  nonzero  byte  is 
found.  This  byte  is  then  scanned  until  the  first  nonzero  bit  is 
found.  The  resulting  exponent  is  computed  as  the  negative  value 
of  eight  times  the  number  of  zero  valued  bytes  plus  the  number  of 
consecutive  zero  valued  bits  in  the  last  byte. 

Expected  Effort 

While  a  fixed  number  of  random  bytes  are  always  used  to  generate 
a  random  mantissa,  the  number  of  byt*  n  used  to  generate  a 
random  exponent  is  itself  a  random  vari  ,  It  can  be  shown  that 


the  distribution  of  n  is: 


PrCn  =  i  >  =  1/2 


for  i  =  0 


i 

=  255  *  (1/256)  /2  for  i  =  1,2,3,... 


Some  of  the  values  of  this  distribution  are  shown  in  Figure  3: 

♦ - - - ♦ - ♦- - ♦ - 1 - ♦ 

!  a  !  0  !  1  !  2  !  3  ! 

♦ - ♦ - ♦ - ♦ - ♦ - ♦ 

!  Prla  bytes  required)  !  .5000  !  .4900408  !  .0019454  !  .0000075  ! 

¥ - ♦ - ♦ - ♦ - ♦ - ♦ 


Figure  3:  Number  of  bytes  required  to  generate  e. 

The  expected  value  of  n  is:  E(n)  =  128  /  255  =  0.5019607 

Assuming  that  one  integer  multiplication  is  required  to  generate 
one  random  byte,  and  that  two  integer  multiplications  are  re¬ 
quired  to  generate  a  random  mantissa,  then  an  average  of  2.5 
integer  multiplications  (  and  no  floating  point  operations)  are 
required  per  floating  point  random  variable.  A  more  efficient 
random  byte  generator  may  bring  this  level  of  effort  down  even 
f urther . 


Expected  period 

The  period  of  the  resulting  generator  is  a  function  of  the 
periods  of  the  generators  used  to  get  the  random  mantissa  m 
(P<m>)  and  the  random  exponent  •  (P(e)>.  If  these  generators 
have  identical  periods,  then  the  random  variable  u  will  also  have 
this  period.  However  if  P<e)  and  P<m>  are  relative  prime  then  the 
period  of  u  i s  equal  to  the  product  of  these  periods. 

P<m)  is  readily  established  as  a  function  of  the  period  of  the 
driving  byte  generator.  However  P(e)  is  more  elusive.  Its  value 
appears  to  depend  on  the  value  of  the  multiplier  used  in  the 
driving  random  byte  generator.  Empirical  tests  using  the  programs 
presented  later  (Figures  7  and  8)  show  that  P(e)  is  practically 
unmeasurable  for  some  values  (3993,  9237  and  14789)  of  a  (the 
multiplier  ,Fig.  8),  while  for  other  values  (1221,2837,2501,7261 
and  12125)  it  is  approx i matel y  equal  to  (but  relative  prime  to) 
the  period  of  the  driving  generator. 

In  both  cases  it  is  reasonable  to  expect  that  the  resulting 
period  of  u  will  be  very  large.  Our  empirical  results  bear  out 
this  conclusion. 


#=»|v|  IMPLEMENTATION 


Notation 

The  AMD  9511  implementation  of  fioatinq  point  notation  is  shown 

in  Figure  4: 

bit  31  30  29  28  25  24  23  14  15  87  0 

♦ - * - - — 1 - +--♦-+ - 

!  sion  !  sian  !  value  of  !  Mantissa  ! 

!  of  !  or  !  of  exponent  least  signif!  I  least  signif! 

! exponent  I eantissa!  !  I  I  I i cant  byte  I  ! leant  byte  I 

♦ - ♦ - - - ♦ - ♦ - ♦ 


Figure  4:  AMD  9511  convention  tor  floating  point  notation 

Three  bytes  are  allocated  to  the  mantissa  and  one  byte  is  used 
for  the  exponent  and  the  two  sign  bits.  Note  that  the  first  bit 
of  the  mantissa  <  always  a  one)  is  explicitly  included  as  bit  24. 
Since  the  first  two  bits  of  byte  1  are  used  for  sign  information, 
only  6  bits  are  available  for  the  order  of  magnitude  o*  the 
exponent.  Exponents  are  therefore  restricted  to  the  range  of  t-64 
to  -*-633. 


For  numbers  between  zero  and  one  the  sign  bit  of  the  mantissa  is 
zero  and  the  sign  bit  of  the  exponent  is  one  if  the  number  is 
less  than  0.5  and  zero  otherwise.  The  effect  of  these  sign  bits 
and  the  fact  that  the  first  bit  of  m  is  always  1  is  shown  in 
Figure  5. 


I  it!  ■  !bl  !  b2  !  b3,b4 

1.0  >*>*0.5  !  0  116,777, 215-8, 388,608  !  0  ’ 

128-2551  0-255 

'0.5  >*>=0.25  1-1  116,777,215-8,388,608  1127' 

128-255!  0-255 

'0.25  >x>=0.125  1-2  U6, 777, 215-8, 388,608  1126’ 

128-255!  0-255 

10.125  >*>=0.0625  !-3  116,777,215-8,388,608  <125 

128-255'  0-255 

!0.0625>x >=0.031251-4  i 16, 777,215-8, 388,608  il24l 

128-255!  0-255 

Figure  5:  Values  of  individual  bytes  for  numbers  on  CO-13 


The  IEEE  implementation  of  this  notation  is  shown  in  Figure  6. 


bit  31  30  29  23  22  21  20  2  1  0 

♦ - ♦ — ♦ — ♦ - - + — , — «. — f 

I  Sion  I  Biased  exponent  I  Unsigned  eantissa  uith  lost  1 
I  of  I  The  bias  is  127  !  significant  bit  tissing  1 

mantissa!  I  !  !  !  !  !  '  !  1  !  ' 

_ i _ x _ x _ x 


It  is  seen  that  the  main  differences  between  this  notation  and 
the  AMD9511  notation  are  the  facts  that  (1)  the  most  significant 
bit  of  the  mantissa  (always  a  one)  is  not  included  in  the  IEEE 
standard  while  it  is  in  the  AMD9511  notation,  (2>  the  IEEE  stan¬ 
dard  allocates  one  more  bit  to  the  exponent,  thus  expanding  the 
range  of  values  that  can  be  represented,  and  (3)  the  exponent- 
mantissa  boundary  is  no  longer  at  a  byte  boundary.  None  of  these 
differences  introduce  serious  implementation  problems.  Note  howe¬ 
ver,  that  the  logical  and  physical  order  of  bytes  in  a  variable 
are  not  always  the  same.  For  example  in  many  systems  the  least 
significant  byte  of  an  integer  is  stored  below  the  more  signifi¬ 
cant  one.  This  may  cause  the  exponent  to  appear  as  the  second 
rather  than  the  first  byte  of  the  variable. 

A  Pascal  Procedure 

The  main  difficulty  in  implementing  the  proposed  algorithm  in 
Pascal  is  to  get  access  to  the  individual  bits  of  a  real  variab¬ 
le.  As  shown  in  the  program  listing  in  Figure  7,  we  obtain  this 
access  through  the  use  of  a  variant  record  that  defines  a  given 
memory  location  to  v  ntain  both  a  four  byte  real  number  and  four 
individual  bytes.  (The  same  effect  is  obtained  in  FORTRAN  through 
the  use  of  an  EQUIVALENCE  statement). 

The  procedure  is  driven  by  the  two  independent  random  byte  gene¬ 
rators  RBytel  and  RByte2.  As  shown  in  Figure  8,  they  obtain  their 
random  bytes  rather  crudely  by  returning  the  first  byte  of  a  two 
byte  integer  produced  by  the  mixed  congruent ial  generator: 

i(k-H)  =  a  *  i<k>  +■  c  mod  (32,768) 

where  :  i <k)  =  i-th  integer  produced  by  the  generator 
c  =  1 

a  =  4  *  j  +  1  <  0  <  j  <  32,768  ) 

It  can  be  shown  C43  that  each  period  of  this  generator  will 
contain  exactly  one  representation  of  all  integers  in  its  range. 
Furthermore,  if  a  is  carefully  selected,  then  the  resulting 
integers  will  pass  most  reasonable  tests  for  randomness  of  dis¬ 
tribution  and  sequence.  Since  these  tests  are  based  on  the  nume¬ 
ric  value  of  the  integer,  it  follows  that  its  most  significant 
byte  of  the  integer  will  also  pass  these  tests.  It  would  be  nice 
if  the  second  byte  also  exhibited  such  properties.  Unf ortunatel  y 
we  were  not  able  to  identify  a  multiplier  a  (we  tried  them  all> 
for  which  this  was  true.  This  is  the  reason  why  the  second  byte 
is  discarded  in  the  generator  shown  in  Figure  8. 

Since  exactly  3  bytes  are  required  to  construct  a  mantissa  while 
an  average  of  0.5019607  bytes  are  required  to  construct  an 
exponent,  the  generators  will  never  be  synchronized.  It  therefore 
seems  quite  safe  to  use  the  same  multiplier  for  both  generators. 
This  was  done  in  the  following  implementation  and  evaluation 
sections. 
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function  unifort  :real; 

typ« 

real type  = 

record  (variant  record  for  byte  access) 

case  integer  of 

1:  (unit  :  real);  (The  variable  ne  want) 

2:  (exponent  :  byte;  (exponent  and  sign  bits) 

al  :  byte;  (tost  significant  byte  of  aantissa) 

•2  :  byte; 

•3  i  byte);  (least  significant  byte  of  tantissa) 
end; 
var 

k  :  byte: 
u  :  real type; 
begin 
uith  u  do 
begin 

■1  :=  RBytel;  (a  separate  byte  generator  is  assuted) 

•2  s*  RBytel; 

•3  :»  RBytel: 

exponent  :=  0:  (  proper  value  for  0.5  <  u  <1.0  > 

if  tl  <  128  then  (by  convention  the  tost  significant  > 

begin  (bit  1  in  al  aust  be  1  > 

al  :*  al  ♦  128: 
exponent  :*  127; 
k  :=  RByte2; 

Nhile  k  *  0  do 
begin 

exponent  ;=  exponent  -  8; 
k  :=  RByte2; 

Hid* 

if  k  J  128  then  begin 
if  k  >=  64  then  exponent  :*  exponent  -1 
else  if  k  >=  32  then  exponent  :=  exponent  *2 
else  if  k  >=  16  then  exponent  :=  exponent  -  3 
else  if  k  >=  8  then  exponent  exponent  -4 
else  if  k  >=  4  then  exponent  :=  exponent  -  5 
else  if  k  >=  2  then  exponent  :-  exponent  -6 
else  exponent  :=  exponent  -  7 
end; 

unrfora  :*  unif; 
end; 
end; 

Figure  7:  Pascal /MT+  Implementation  of  uniform  random  number 
generator  using  AMD  9511  floating  point  notation. 


function  RBytel:Byte; 

(  note:  the  folloaing  declarations  aust  appear  in  the  aain  prograa 
VAR 

seed  :  record 

case  integer  of 
1:  (int  :  integer); 

2:  <lsbyte:byte: 

asbyte:byte); 
end;  } 

CONST 

MULT  -  1221;  (  other  good  values  are:  2837,  3993,  4184,  4293.) 
begin  (  9237,14789,15123,17245) 

seed.int  :=  WILT  t  seed. int  ♦  *; 

RandoaByte  :=  seed.asbyte; 
if  seed.int  <  0  then  seed.int  *»nt*l; 

end; 


Figure  8.  Random  byte  generator  for  Z80  based  system.  Note 

that  the  order  of  bytes  in  the  integer  is  non  standard 


I  ON 

A  summary  of  the  tests  performed  to  assess  the  quality  of  the 
proposed  generator  is  presented  in  Figure  9.  For  each  test  the 
initial  seed  tor  RBytel  was  12345  and  the  initial  seed  for  RByte2 
was  2345.  To  save  space  details  of  the  test  procedures  are  not 
presented  here.  The  reader  is  referred  to  C&3  or  to  the  reference 
cited  in  Figure  9  for  more  information. 


!  Property  !  Eapirical  Measure 

Criteria  !  Procedure  !  Reference  ! 

!  Unifora  !  100  bin  frequency  count 

!  Distribution  ! 

Chi-square  !  1000  observations  !  14]  p.35  ! 

99  df  !  1000  replications  !  ! 

!  Randoa  !  Count  of  runs  up  of  length 

!  Sequence  !  1,2, 3, 4  and  5  or  aore 

Chi-square  !  1000  observations  !  141  p.68  ! 

5  df  !  1000  replications  !  ! 

!  Lack  of  auto  !  Counts  of  noraalized  acf’s 
!  correlation  !  outside  ♦ 1-2  s  for  all  lags 

No  clear  pattern!  2000  observations  !  161  ! 

Most  counts  (4  !  20  replications  !  ! 

!  Coaputational  !  Tiae  to  generate  10,000 
!  Effort  !  variables 

Other  !  arist  watch  tiaing!  ! 

generators  !  on  4Mz  Z80  aicro  !  ! 

!  Period  !  Bap  between  two  identical 

!  !  sequences  of  100  nuabers 

Other  !  Published  results  !  ! 

Generators  !  used  if  available  !  ! 

Figure  9.  Summary  of  Evaluation  tests. 


The  outcome  of  these  when  tests  applied  to  the  generator  in 
Figure  7  using  the  byte  generator  in  Figure  8  are  shown  in  Figure  10 


1221 


HULT 


DISTRIBUTION 


Range  of  Chi-square 
values,  11000  tests) 


■in  !  avg  !  aax 


:==s==fsr=ss==+==s==== 

61.0  !  99.9  .»  155.4 


;+=s==s===t=======5+:====s=4.==s=ssf: 


Test  outcoaes 
at  51  level 


accept 

932 


reject 
68 


RUN-TEST 


nwif  nv  ■  nvininohii  ■ 


Range  of  Chi-square!  Test  outcoaes 
values  11000  obs)  !  at  51  level 


Index  of  lags  w/  too  aany! 
acf’s  outside  ♦ /-  2  sigaa! 


am 

0.1 


avg 

4.i 


■ax  !  accept  !  reject 
20.9  !  979  !  21 


AUTOCORRELATION 


4  bad  acf’s 
109 


5  bad  acf’s! 
none  ! 


— ♦ 


2837 


59.2  !  96.4  !  138.0 


974 


26 


0.1 


4.1 


25.7  !  952  !  42 


51 


none 


3993 


62.6  !  99.0  !  163.2 


957 


43 


0.1 


3.8 


19.1 


979 


21 


none 


4189 


63.0  )  96.0  !  141.0 


976 


24 


0.0 


4.1 


18.1 


983 


17 


none 


17 


4293 


57.2  1  95.0  !  138.0 


976 


24 


0.0 


3.8 


20.5  !  985  !  15 


41 


77 


9237 


66.4  !  99.6  !  151.2 


936 


64 


0.1 


4.0 


23.1  !  971  !  29 


61 


none 


14789 


61.2  !  98.6  !  150.4 


958 


42 


0.1 


4.2 


17.9  !  971 


29 


20,  65 


none 


15125 


69.0  !  98.5  !  146.8 


979 


21 


0.1 


4.0 


14,9 


977 


23 


1,  115 


none 


!  17245 


67.2  '  99.5  !  137.2 


!  IDEAL 


66.5  !  98.3  !  139.0 


960 
50 


40 

50 


0.1 


4.0 


0.4 


4.3 


19.4  !  974  !  26 
16.6  !  950  !  50 


25 


2  bad  lags 


42  ! 

erssssssssss s+ 

no  bad  lags  1 


Figure  10s  Results  of  Empirical  tests  for  recommended  seeds. 
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It  should  be  noted  that  our  subjective  criteria  of  "closeness  tc 
Chi-square  distribution"  is  a  more  discriminating  test  than  the 
conventional  chi-square  test  discussed  in  L4I  .  While  the  convert 
tional  test  rejects  the  hypothesis  it  a  single  value  ot  the  test 
statistic  is  out  ot  bound,  we  generate  many  values  o+  the  statis¬ 
tic  and  reject  if  the  resulting  empirical  distribution  does  not 
appear  to  follow  the  Chi-square  dis  tribution.  (We  dre  able  to  do 
this  since  the  long  period  of  the  generator  allows  us  access  to 
an  extremely  large  samplesize) 


In  order  to  determine  the  effort  required  to  generate  a  t  -tndom 
number,  we  measured  the  time  required  to  generate  10, 000 
different  numbers  using  the  present  and  other  algorithms.  All 
programs  were  run  on  a  4Mz  Z80  based  micro  computer  without 
special  floating  point  hardware.  The  results  ot  these  tests  are 
summarized  in  Figure  11; 


Setter  at  or  type 

l 

!  period 

resolution  !n 

ip  between  1  1 

ke  numbers'  Tite  ' 

Figure  7 

1  T  V 

!  at  least  !  -10! 

164,000,000!  4.7  1  10  1 

randoi  !  18  ! 

1  1 

Figure  7  using  a  Tausuortne  type 

!  at  least 

-10' 

randoi  14  1 

randoi  byte  generator  13,5,61 


164,000,000!  4.7  »  10 


Multiplicative  congruential  generator  1 
returning  a  16  bit  integer  !  32,768 

I(n+1)  =  1  In)  *  a  I  1  Hod  32768  1 


32,768 


!  Multiplicative  congruential  generator 
'  returning  a  value  on  0-1 


!  U  =  I <n) 10. 00003051757 


32,3278 


-5 


3.1  *  10 


32,786  !  18 


!  Multiplicative  congruential  generator  ! 

!  emulating  32  bit  arithmetic  in  16  bits!  32,3278 
!  Ilntb  =  1  <n>  I  a  *  1  Mod  21*31  ! 

♦ - - - 


32,786  1  133 


Figure  11.  Comparison  between  different  generators. 


Two  implementations  of  the  prosed  generator  were  tinned.  The 
first  [Figure  73  used  the  random  byte  generator  presented  in 
Figure  8,  while  the  second  used  Figure  7  in  conjunction  with  an 
experimental  implementation  of  a  Tausworthe  type  random  byte 
generator  13,5,61.  The  first  of  these  was  about  as  fast  as  a 
conventional  LCG  based  floating  point  generator  generator,  while 
the  second  was  substantially  faster.  Both  generators  yielded 
random  number  streams  with  statistical  properties  vastly  superior 
to  LCG  based  generators. 
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SUMMARY 


A  procedure  tor  the  construction  ot  uniformly  distributed 
deviates  on  CO-ID  -from  two  independently  distributed  streams  of 
random  numbers  has  been  presented.  Based  on  the  conventions  used 
■for  floating  point  notation,  the  procedure  is  shown  to  have  good 
statistical  properties.  This  includes  the  fact  that  the  interval 
between  two  like  numbers  is  itself  a  random  variable  (the  absence 
of  this  property  is  frequently  a  reason  for  concern  when  conven¬ 
tional  congruential  generators  are  used).  The  procedure  is  shown 
to  be  computationally  efficient,  and  further  improvements  in  com¬ 
putational  speeds  should  be  possible  through  the  use  of  a  more 
efficient  generator  of  random  bytes.  One  disadvantage  of  the 
procedure  is  a  lack  of  portability  of  the  resulting  computer 
proqrams.  This  is  because  (1)  different  machines  used  different 
conventions  for  representing  floating  point  numbers,  and  (2)  the 
driving  byte  generators  are  likely  to  be  machine  dependent 
congruential  generators.  However  these  restrictions  are  about  the 
same  as  those  experienced  for  conventional  congruential 
generator s. 
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